function Yest=rbf_eval(X,S,lambda,gamma,flag)
% 用于RBF验证
%
% INPUTS:
%   X      - 需要预测函数值的点 (待估计点)
%   S      - 已知函数值的采样点 (训练样本)
%   lambda - 通过拟合得到的权重参数向量
%   gamma  - 多项式尾项的参数向量
%   flag   - 指定所使用的径向基函数类型:
%           'cubic'        : 三次径向基函数 φ(r) = r^3
%           'TPS'          : 薄板样条 φ(r) = r^2 * log(r)
%           'linear'       : 线性径向基函数 φ(r) = r
%           'Gaussian'     : 高斯径向基函数 φ(r) = exp(-(r/ρ)^2)
%           'multiquad'    : 多二次径向基函数 φ(r) = sqrt(r^2 + ρ^2)
%           'invmultiquad' : 反多二次径向基函数 φ(r) = 1/sqrt(r^2 + ρ^2)
%
% OUTPUT:
%   Yest   - 在 X 上预测的函数值

[mX,nX]=size(X); % 带预测点的维度
[mS,nS]=size(S); % 已知采样点的维度
if nX~=nS % 检查维度是否一致，如果不一致则转置 X
    X=X';
    [mX,nX]=size(X);
end

R=zeros(mX,mS); % 存储 X 与 S 的两两距离
% 优先使用 MATLAB 自带的 pdist2 计算欧氏距离
if exist('pdist2', 'file')
    R = pdist2(X,S,'euclidean');
else
    % 如果没有 pdist2，则手动计算欧氏距离
    for i = 1:mX
        for j = 1:mS
            R(i,j) = sqrt(sum((X(i,:) - S(j,:)).^2));
        end
    end
end

% 根据 flag 选择不同的径向基函数
if strcmp(flag,'cubic')
    Phi= R.^3;
elseif strcmp(flag,'TPS')
    R(R==0)=1; % 避免 log(0)
    Phi=R.^2.*log(R);
elseif strcmp(flag,'linear')
    Phi=R;
elseif strcmp(flag,'Gaussian')
    rho = max(max(R))/(nX*mS)^(1/nX); % 自适应参数 ρ
    Phi=exp(-R.^2/rho^2);
elseif strcmp(flag,'multiquad')
    rho = max(max(R))/(nX*mS)^(1/nX);
    Phi=sqrt(R.^2 + rho^2);
elseif strcmp(flag,'invmultiquad')
    rho = max(max(R))/(nX*mS)^(1/nX);
    Phi=1./sqrt(R.^2 + rho^2);
end

% 预测值由两部分组成：RBF 主体 + 多项式尾项
Yest1=Phi*lambda;          % RBF 主体
Yest2=[X,ones(mX,1)]*gamma; % 多项式尾项
Yest=Yest1+Yest2;           % 最终预测函数值
end % function
